We removed “one-off” species as defined by appearing less than five times across all species datasets and if the parasite also exists in its own group (e.g., We found four acanthacephalean species across four fish species, so acanths were removed because less than five were found and they exist in a seperate group from other parasites. e.g., TREM.META was kept because even though it only appeared once in the H. nuchalis dataset, it will be lumped in with the other trematode metacercaria data points.)
| Variable | Definition |
|---|---|
| K | Fulton’s Condition Factor (<1=good, >1=bad) |
| W | weight of host fish (in grams) |
| SL | standard length of fish (in millimeters) |
\[ K = W/SL^3 *100 \]
library(ggeffects)
library(tidyverse)
library(ggpubr)
library(caret)
library(tibble)
library(readr)
library(dbplyr)
library(lubridate)
library(ggplot2)
library(glmmTMB)
library(DHARMa)
H.nuchalis.data <- H.nuchalis.data %>%
mutate(weight_mg = as.numeric(weight_mg))
Pvig.life.stage.total <- P.vigilax.data %>%
mutate(adult = rowSums(select(.,
"acanth_bigb",
"acanth_bkr",
"nem_cona",
"nem_cap",
"nem_dich",
"nem_larv",
"nem_trbk",
"nem_unk",
"cest_godz",
"cest_comp",
"cest_ntg",
"cope_cl",
"cope_grab",
"cest_pea",
"cest_l",
"cest_vit",
"cest_unk",
"mono_gyro",
"mono_dact",
"mono_lg",
"trem_rnda",
"trem_unk"))) %>%
mutate(larval = rowSums(select(.,
"meta_hs",
"meta_unk",
"myx_go",
"myx_geom",
"myxo_bc",
"myx_thel",
"myxo_sbad",
"myxo_unk",
"myxo_up",
"trem_met",
"trem_sss",
"trem_buc",
"trem_metaf",
"trem_metags",
"trem_ms",
"trem_nea"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
janitor::clean_names()
####
Ipun.life.stage.total <- I.punctatus.data %>%
mutate(larval = rowSums(select(.,
"nem_cyst",
"trem_meta_gg",
"trem_meta_unk",
"trem_lg",
"nem_cyst",
"myx_tail",
"myx_geom"))) %>%
mutate(adult = rowSums(select(.,
"acan_ostr",
"cest_comp",
"cest_meg",
"leech_kut",
"nem_nmts",
"mono_ip",
"mono_unk",
"nem_phar",
"nem_pty",
"nem_sp",
"nem_taf",
"nem_unk",
"nmorph",
"trem_2",
"trem_allo",
"trem_ble",
"trem_crep",
"trem_f",
"trem_lip",
"trem_megict",
"trem_smile",
"trem_unk"))) %>%
select(-nmorph, -leech_kut) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
janitor::clean_names()
####
Nath.life.stage.total <- N.atherinoides.data %>%
mutate(larval = rowSums(select(.,
"nem_cyst",
"trem_larv",
"trem_meta",
"myx_sp",
"myx_round",
"trem_cyst"))) %>%#
mutate(adult = rowSums(select(.,
"cest_comp",
"cest_hyd",
"cest_lvs",
"cest_tri",
"cope_a",
"mono_all",
"mono_unk",
"nem_bran",
"nem_cap",
"nem_spky",
"nem_tfs",
"nem_twas",
"nem_unkn",
"trem_gb",
"trem_l"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
janitor::clean_names()
####
Hnuc.life.stage.total <- H.nuchalis.data %>%
mutate(larval = rowSums(select(.,
"nem_daft",
"trem_cm",
"trem_meta_het",
"trem_meta_sp",
"trem_bc",
"trem_meta_go",
"trem_meta_es",
"meta_unk",
"trem_meta",
"trem_meta_unk",
"myx_ak",
"myx_fi",
"myx_tin",
"myx_gl",
"myx_ot",
"myx_eye",
"myx_gc"))) %>%#
mutate(adult = rowSums(select(.,
"cest_bb",
"nem_bud",
"nem_unk",
"trem_ass",
"trem_ns",
"trem_diplo",
"trem_unk",
"cest_unk",
"mono_dact",
"mono_gdac",
"acan_ad",
"cope_pood"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
janitor::clean_names()
####
Pvig.life.stage.pared <- P.vigilax.data %>%
mutate(adult = rowSums(select(.,
"acanth_bigb",
"acanth_bkr",
"nem_cona",
"nem_cap",
"nem_dich",
"nem_larv",
"nem_trbk",
"nem_unk",
"cest_godz",
"cest_comp",
"cest_ntg",
"cope_cl",
"cope_grab",
"cest_pea",
"cest_l",
"cest_vit",
"cest_unk",
"mono_gyro",
"mono_dact",
"mono_lg",
"trem_rnda",
"trem_unk"))) %>%
mutate(larval = rowSums(select(.,
"meta_hs",
"meta_unk",
"myx_go",
"myx_geom",
"myxo_bc",
"myx_thel",
"myxo_sbad",
"myxo_unk",
"myxo_up",
"trem_met",
"trem_sss",
"trem_buc",
"trem_metaf",
"trem_metags",
"trem_ms",
"trem_nea"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
select(-matches("^(myx_|myxo_|acanth_|mono_|cope_)")) %>%
janitor::clean_names()
####
Ipun.life.stage.pared <- I.punctatus.data %>%
mutate(larval = rowSums(select(.,
"nem_cyst",
"trem_meta_gg",
"trem_meta_unk",
"trem_lg",
"nem_cyst",
"myx_tail",
"myx_geom"))) %>%
mutate(adult = rowSums(select(.,
"acan_ostr",
"cest_comp",
"cest_meg",
"leech_kut",
"nem_nmts",
"mono_ip",
"mono_unk",
"nem_phar",
"nem_pty",
"nem_sp",
"nem_taf",
"nem_unk",
"nmorph",
"trem_2",
"trem_allo",
"trem_ble",
"trem_crep",
"trem_f",
"trem_lip",
"trem_megict",
"trem_smile",
"trem_unk"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
select(-matches("^(myx_|acan_|mono_|leech_|nmorph)")) %>%
janitor::clean_names()
####
Nath.life.stage.pared <- N.atherinoides.data %>%
mutate(larval = rowSums(select(.,
"nem_cyst",
"trem_larv",
"trem_meta",
"myx_sp",
"myx_round",
"trem_cyst"))) %>%#
mutate(adult = rowSums(select(.,
"cest_comp",
"cest_hyd",
"cest_lvs",
"cest_tri",
"cope_a",
"mono_all",
"mono_unk",
"nem_bran",
"nem_cap",
"nem_spky",
"nem_tfs",
"nem_twas",
"nem_unkn",
"trem_gb",
"trem_l"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
select(-matches("^(myx_|mono_|cope_)")) %>%
janitor::clean_names()
####
Hnuc.life.stage.pared <- H.nuchalis.data %>%
mutate(larval = rowSums(select(.,
"nem_daft",
"trem_cm",
"trem_meta_het",
"trem_meta_sp",
"trem_bc",
"trem_meta_go",
"trem_meta_es",
"meta_unk",
"trem_meta",
"trem_meta_unk",
"myx_ak",
"myx_fi",
"myx_tin",
"myx_gl",
"myx_ot",
"myx_eye",
"myx_gc"))) %>%#
mutate(adult = rowSums(select(.,
"cest_bb",
"nem_bud",
"nem_unk",
"trem_ass",
"trem_ns",
"trem_diplo",
"trem_unk",
"cest_unk",
"mono_dact",
"mono_gdac",
"acan_ad",
"cope_pood"))) %>%
mutate(adult = as.integer(adult)) %>%
mutate(larval = as.integer(larval)) %>%
mutate(psite_count = rowSums(select(.,
"adult",
"larval"))) %>%
select(-matches("^(myx_|mono_|cope_|acan_)")) %>%
janitor::clean_names()
### For answering question 1
Nath.life.stage.total <- Nath.life.stage.total %>%
mutate(body.index.ath = ((Nath.life.stage.total$weight_mg/((Nath.life.stage.total$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
Pvig.life.stage.total <- Pvig.life.stage.total %>%
mutate(body.index.vig = ((Pvig.life.stage.total$weight_mg/((Pvig.life.stage.total$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
Ipun.life.stage.total <- Ipun.life.stage.total %>%
mutate(body.index.pun = ((Ipun.life.stage.total$weight_mg/((Ipun.life.stage.total$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
Hnuc.life.stage.total <- Hnuc.life.stage.total %>%
mutate(body.index.nuc = ((Hnuc.life.stage.total$weight_mg/((Hnuc.life.stage.total$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
### For answering question 2
Nath.life.stage.pared <- Nath.life.stage.pared %>%
mutate(body.index.ath = ((Nath.life.stage.pared$weight_mg/((Nath.life.stage.pared$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
Pvig.life.stage.pared <- Pvig.life.stage.pared %>%
mutate(body.index.vig = ((Pvig.life.stage.pared$weight_mg/((Pvig.life.stage.pared$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
Ipun.life.stage.pared <- Ipun.life.stage.pared %>%
mutate(body.index.pun = ((Ipun.life.stage.pared$weight_mg/((Ipun.life.stage.pared$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
Hnuc.life.stage.pared <- Hnuc.life.stage.pared %>%
mutate(body.index.nuc = ((Hnuc.life.stage.pared$weight_mg/((Hnuc.life.stage.pared$standard_length_mm)^3)))*(100)) %>%
janitor::clean_names()
# plot condition factor for life stage (aka "stage-specific parasite stuff")
NA.bc.a <- ggplot(Nath.life.stage.pared,
aes(x = adult, y = body_index_ath, color = sex)) + xlab("Adult") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
NA.bc.l <- ggplot(Nath.life.stage.pared,
aes(x = larval, y = body_index_ath, color = sex)) + xlab("Larval") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
# plot condition factor for life stage (aka "stage-specific parasite stuff")
IP.bc.a <- ggplot(Ipun.life.stage.pared,
aes(x = adult, y = body_index_pun, color = sex)) + xlab("Adult") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
IP.bc.l <- ggplot(Ipun.life.stage.pared,
aes(x = larval, y = body_index_pun, color = sex)) + xlab("Larval") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
# plot condition factor for life stage (aka "stage-specific parasite stuff")
PV.bc.a <- ggplot(Pvig.life.stage.pared,
aes(x = adult, y = body_index_vig, color = sex)) + xlab("Adult") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
PV.bc.l <- ggplot(Pvig.life.stage.pared,
aes(x = larval, y = body_index_vig, color = sex)) + xlab("Larval") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
# plot condition factor for life stage (aka "stage-specific parasite stuff")
HN.bc.a <- ggplot(Hnuc.life.stage.pared,
aes(x = adult, y = body_index_nuc, color = sex)) + xlab("Adult") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
HN.bc.l <- ggplot(Hnuc.life.stage.pared,
aes(x = larval, y = body_index_nuc, color = sex)) + xlab("Larval") + geom_point() +
ylab("Fulton's Condition Factor") + theme_bw() + stat_smooth(method=lm)
ggarrange(IP.bc.a, IP.bc.l, NA.bc.a, NA.bc.l, PV.bc.a, PV.bc.a, HN.bc.a, HN.bc.l, common.legend = TRUE, legend = "bottom", labels = c("A", "B", "C", "D", "E", "F", "G", "H"))
# more variability than the others
## double check sex, weights
## should we exclude these outliers
We want to run a general linear model to test the research questions: 1. Is body condition associated with increased total parasite burden? - is there a linear association between increasing body condition and increasing parasite burden? - does FCF increase as parasite burden increases? 2. Is body condition associated with infection with adult or larval stages of parasites? - is there a linear association between increasing body condition and increasing abundance of larval or adult parasites in the host? - does FCF increase as adult or larval parasite presence increases? - variables: - predictor = parasite_count, life_stage –> discrete - response = FCF –> continuous - discrete predictor = ANOVA - multiple/ continuous and discrete predictors = ANCOVA?
\[ K_i\sim parasite\ count \]
\[ K_i\sim\ life\ stage\ + parasite\ count \] #### Model formula for Question 3 \[ K_i\sim\ parasite\ stage\ * parasite\ group\ + parasite\ count \]
# Question 1 (lm)
summary(FCF_ath_lm <- lm(body_index_ath ~ psite_count, data = Nath.life.stage.total))
##
## Call:
## lm(formula = body_index_ath ~ psite_count, data = Nath.life.stage.total)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.8708 -0.2325 -0.1553 -0.0608 7.2784
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.088027 0.085561 12.716 <2e-16 ***
## psite_count 0.001006 0.003980 0.253 0.801
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9823 on 169 degrees of freedom
## (37 observations deleted due to missingness)
## Multiple R-squared: 0.0003782, Adjusted R-squared: -0.005537
## F-statistic: 0.06395 on 1 and 169 DF, p-value: 0.8007
plot(FCF_ath_lm)
summary(FCF_ath_lm <- lm(body_index_pun ~ psite_count, data = Ipun.life.stage.total))
##
## Call:
## lm(formula = body_index_pun ~ psite_count, data = Ipun.life.stage.total)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.8625 -0.5370 -0.3616 -0.1666 13.7487
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.860079 0.261188 7.122 6.59e-10 ***
## psite_count -0.002641 0.005931 -0.445 0.657
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.14 on 72 degrees of freedom
## (15 observations deleted due to missingness)
## Multiple R-squared: 0.002746, Adjusted R-squared: -0.01111
## F-statistic: 0.1982 on 1 and 72 DF, p-value: 0.6575
plot(FCF_ath_lm)
summary(FCF_ath_lm <- lm(body_index_vig ~ psite_count, data = Pvig.life.stage.total))
##
## Call:
## lm(formula = body_index_vig ~ psite_count, data = Pvig.life.stage.total)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.8966 -0.4336 -0.2843 -0.0619 12.7554
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.997946 0.305819 6.533 3.51e-08 ***
## psite_count -0.003891 0.008464 -0.460 0.648
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.874 on 49 degrees of freedom
## (146 observations deleted due to missingness)
## Multiple R-squared: 0.004294, Adjusted R-squared: -0.01603
## F-statistic: 0.2113 on 1 and 49 DF, p-value: 0.6478
plot(FCF_ath_lm)
summary(FCF_nuc_lm <- lm(body_index_nuc ~ psite_count, data = Hnuc.life.stage.total))
##
## Call:
## lm(formula = body_index_nuc ~ psite_count, data = Hnuc.life.stage.total)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.46132 -0.14780 0.01483 0.14594 1.36842
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.613e+00 2.419e-02 66.657 <2e-16 ***
## psite_count -7.604e-05 2.325e-04 -0.327 0.744
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2874 on 201 degrees of freedom
## (18 observations deleted due to missingness)
## Multiple R-squared: 0.0005318, Adjusted R-squared: -0.004441
## F-statistic: 0.107 on 1 and 201 DF, p-value: 0.744
plot(FCF_nuc_lm)
# Question 2 (lm)
summary(FCF_ath_al_lm <- lm(body_index_ath ~ larval + psite_count, data = Nath.life.stage.pared))
##
## Call:
## lm(formula = body_index_ath ~ larval + psite_count, data = Nath.life.stage.pared)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.9131 -0.2363 -0.1636 -0.0531 7.2716
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.076895 0.092623 11.627 <2e-16 ***
## larval -0.003248 0.010187 -0.319 0.750
## psite_count 0.003675 0.009273 0.396 0.692
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9849 on 168 degrees of freedom
## (37 observations deleted due to missingness)
## Multiple R-squared: 0.0009826, Adjusted R-squared: -0.01091
## F-statistic: 0.08262 on 2 and 168 DF, p-value: 0.9207
plot(FCF_ath_al_lm)
summary(FCF_ath_al_lm <- lm(body_index_pun ~ larval + psite_count, data = Ipun.life.stage.pared))
##
## Call:
## lm(formula = body_index_pun ~ larval + psite_count, data = Ipun.life.stage.pared)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.8958 -0.5067 -0.3769 -0.1667 13.6968
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.911985 0.295691 6.466 1.1e-08 ***
## larval 0.009953 0.026007 0.383 0.703
## psite_count -0.011973 0.025104 -0.477 0.635
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.152 on 71 degrees of freedom
## (15 observations deleted due to missingness)
## Multiple R-squared: 0.004799, Adjusted R-squared: -0.02324
## F-statistic: 0.1712 on 2 and 71 DF, p-value: 0.843
plot(FCF_ath_al_lm)
summary(FCF_ath_al_lm <- lm(body_index_vig ~ larval + psite_count, data = Pvig.life.stage.pared))
##
## Call:
## lm(formula = body_index_vig ~ larval + psite_count, data = Pvig.life.stage.pared)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.9023 -0.4604 -0.2923 -0.0385 12.7392
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.03633 0.34361 5.926 3.25e-07 ***
## larval 0.02669 0.10483 0.255 0.800
## psite_count -0.02969 0.10168 -0.292 0.772
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.892 on 48 degrees of freedom
## (146 observations deleted due to missingness)
## Multiple R-squared: 0.005637, Adjusted R-squared: -0.03579
## F-statistic: 0.1361 on 2 and 48 DF, p-value: 0.8731
plot(FCF_ath_al_lm)
summary(FCF_nuc_al_lm <- lm(body_index_nuc ~ larval + psite_count, data = Hnuc.life.stage.pared))
##
## Call:
## lm(formula = body_index_nuc ~ larval + psite_count, data = Hnuc.life.stage.pared)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.45556 -0.14240 0.00991 0.14511 1.36850
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.6198938 0.0263203 61.545 <2e-16 ***
## larval 0.0009821 0.0013692 0.717 0.474
## psite_count -0.0010150 0.0013297 -0.763 0.446
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2877 on 200 degrees of freedom
## (18 observations deleted due to missingness)
## Multiple R-squared: 0.003096, Adjusted R-squared: -0.006873
## F-statistic: 0.3106 on 2 and 200 DF, p-value: 0.7334
plot(FCF_nuc_al_lm)
#df <- df %>%
# group_by(fish_ID, psite_group, life_stage) %>%
#sum(???)